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ABSTRACT 

We use numerical simulations of cosmic reionization and radiative processes related to 
the HI 21 cm emission line to produce synthetic radio maps as seen by next generation 
telescopes that will operate at low radio frequencies (e.g. LOFAR). Two different 
scenarios, in which the end of reionization occurs early {z w 13) or late {z ^ 8) 
depending on the Initial Mass Function (IMF) of the first stars and ionizing photon 
escape fraction, have been explored. For each of these models we produce synthetic 
HI 21 cm emission maps by convolving the simulation outputs with the provisional 
LOFAR sampling function in the frequency range 76-140 MHz. If reionization occurs 
late, LOFAR will be able to detect individual HI structures on arcmin scales, emitting 
at a brightness temperature of w 35 mK as a Z-a signal in about 1000 hours of 
observing time. In the case of early reionization, the detection would be unlikely, due 
to decreased sensitivity and increased sky temperatures. These results assume that 
ionospheric, interference and foreground issues are fully under control. 
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1 INTRODUCTION 



Cosmic recombination left the gas in the universe in a nearly 
uniform, cold and neutral state. This vast sea of baryonic 
material was only rippled by tiny density fluctuations which, 
under the action of gravity, later grew into bound objects in- 
side which eventually stars formed. This evolutionary trend 
was however counteracted by a number of physical processes: 
the finite, albeit low, Jeans mass scale (M « 10* M©) pre- 
venting the gas to closely track the collapse of the underly- 
ing dark matter component; the lack of the standard cooling 
agents (heavy elements, dust) for the gas; the fragility and 
low cooling efficiency of available molecules (H2, HD, LiH); 
finally, a number of radiative and mechanical feedback ef- 
fects (e.g. Ciardi & Ferrara 2005) arising from the formation 
of the first luminous sources. As a result of these complica- 
tions, shaping the universe as we currently monitor it up to 
the most remote objects was a somewhat grueling process. 
Therefore, a non-negligible fraction (« 1%) of cosmic time 
was spent in a relatively dark and quiescent state that, fol- 
lowing M. Rees suggestion, we currently denote as the Dark 
Ages of the universe. 

The investigation of the Dark Ages is still in its infancy, 
and very little, if any, experimental support to theory is cur- 
rently available. However, the most promising (by far) tech- 
nique to explore the evolution of the intergalactic medium 
(IGM) during these remote times is the study of the 21 cm 



hyperfine triplet-singlet level transition of the ground state 
of neutral hydrogen. This line could, in principle, allow a 
superb tracing of the HI distribution in the early universe, 
and therefore a reconstruction of the reionization history as 
governed by the first luminous sources. 

In recent years two independent sets of data coming 
from QSO absorption line experiments, and the analysis 
of the Cosmic Microwave Background (CMB) temperature- 
polarization cross-correlation spectrum, constrained the 
epoch of reionization (EoR). In brief, the sharp raise of the 
Gunn-Peterson Lya (and Ly/3) optical depth derived from 
the analysis of z 5 quasar spectra, has led some authors 
to interpret this result as an indication of the rapid neu- 
tral/ionized transition to take place at relatively low redshift 
(e.g. White et al. 2003). The first year Wilkinson Microwave 
Anisotropy Probe ( WMAP) data instead seemed to indicate 
that either the EoR was located at a higher redshift or that 
a more complex evolution of the gas ionized fraction took 
place (e.g. Kogut et al. 2003). This interpretation has to be 
revised in view of the newly released data (Spergel et al. 
2006), which favor a lower reionization redshift. In recent 
years a large effort has been made to numerically simulate 
the formation of the first structures and the resulting pho- 
toionization of the IGM (e.g. Gnedin 2000; Ciardi, Ferrara, 
& White 2003 (CFW); Sokasian et al. 2003; Gnedin 2004; 
Iliev et al. 2006). Although it is possible to devise models 
which are in agreement with both sets of data (e.g. Choud- 
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hury & Ferrara 2005; Gallerani, Choudhury & Ferrara 2005), 
independent measurements are eagerly required. 

Redshiftcd 21 cm line emission from the Dark Ages rep- 
resents a unique chance to fully map the spatial distribution 
of intergalactic hydrogen (e.g. Madau, Meiksin & Rcos 1997; 
Ciardi & Madau 2003 (CM); Furlanetto, Sokasian & Hern- 
quist 2004; Furlanetto, ZaJdarriaga & Hernquist 2004) as a 
function of redshift. The size of the structures that could be 
detected depends on the design of future radio telescopes. 
The next generation of radio telescopes such as the Square 
Kilometer Array (SKA), the LOw Frequency ARray (LO- 
FAR), the PrimevAl Structure Telescope (PAST) and the 
Mileura Wide-field Array (MWA) will probably have the 
sensitivity for HI mapping at resolutions of the order of a 
few arcminutes (e.g. Pen, Wu, & Peterson 2004; Bowman, 
Morales, & Hewitt 2005; Kassim et al. 2004; Wyithe, Loeb, 
& Barnes 2005). Our aim here is limited, in the sense that 
we will describe the likelihood of one particular instrument 
(i.e. LOFAR) to be able to detect a signal from EoR. To 
do this we have simulated observations of the EoR signal 
for two very general (early, late) reionization scenarios with 
an idealized LOFAR array, using state-of-art cosmological 
radiative transfer numerical computations. Although this is 
only a first step towards a proper modelling, the results al- 
ready give a clear fiavor of the tremendous potential that 
these types of experiments will provide in the very near fu- 
ture. 



The 21 cm radiation intensity can be expressed by the 
differential brightness temperature between a neutral hydro- 
gen patch and the CMB: 



5Tt 



Ts — TcML 

1 + z 



(3) 



where r is the optical depth of the neutral IGM at 21(1 -|- 
z) cm: 



_ Sc^hpAio 



(4) 



In eq. 4, hp and k are the Planck and Boltzmann constants, 
respectively; vq = 1420 MHz is the 21 cm hyperfine transi- 
tion frequency. Aw = 2.85x 10~^^s~^ is the 21 cm transition 
spontaneous decay rate, and nui is the local HI number den- 
sity. 

If Tk (and thus Ts) is higher than Tcmb the neutral 
IGM will be visible in emission against the CMB; on the 
contrary, if Ts < Tcmb it will be visible in absorption. The 
heating by primordial sources of radiation should easily push 
Ts to values ^ Tcmb at the redshifts of interest (see e.g. 
Madau, Meiksin & Rees 1997; Carilli et al. 2002; Venkatesan, 
Giroux & ShuU 2001; Chen & Miralda-Escude 2004). Then 
the IGM will be observable in emission at a level that is 
independent of the exact value of Ts and is directly propor- 
tional to riHi- A radio interferometer able to make a 21 cm 
tomography would allow to probe accurately the reionizar 
tion history. 



2 21 CM RADIO EMISSION 

The 21 cm lino is associated with the hyperfine transition 
between the triplet and the singlet levels of the hydrogen 
ground state. This transition is governed by the spin tem- 
perature, Ts, defined as: 




where no and m are the number densities of hydrogen atoms 
in the singlet and triplet ground hyperfine levels, and T^, = 
0.068 K is the temperature corresponding to the transition 
energy. 

In the presence of the CMB radiation alone the spin 
temperature reaches thermal equilibrium with Tcmb = 
2.73(1 -I- z) K on a short timescale, making the HI unde- 
tectable in emission or absorption. A mechanism is then re- 
quired to decouple Ts from Tcmb to observe the line. For 
the range of densities typical of the IGM, the scattering by 
Lya photons is an efficient mechanism, as it couples Ts to 
the kinetic gas temperature Tk mixing the hyperfine levels 
of the neutral hydrogen in its ground state through interme- 
diate transitions to the excited 2p state. This is the so called 
Wouthuysen-Field process, or Lya pumping (e.g. Wouthuy- 
sen 1952; Field 1959; Hirata 2005). For the Lya pumping to 
be effective the Lya background intensity, Ja , has to satisfy 
the condition (e.g. CM): 

Jq > 9 X 10"^-''(1 + z) ergs cm"^ s~^ Hz"^ sr"\ (2) 

at the redshift of interest. CM find that this condition on 
the diffuse Lya intensity is satisfied between 2; w 20 and the 
EoR. 



3 SIMULATIONS OF 21 CM LINE EMISSION 

In this Section we briefly describe the simulations of 21 cm 
line emission devised to produce synthetic brightness tem- 
perature maps; we refer the reader to the original papers 
(Ciardi, Stoehr & White 2003; CFW; CM) for more details. 
The simulations of reionization* employ a combination of 
high resolution N-body simulations (to describe the distri- 
bution of dark matter and diffuse gas; Springel, Yoshida, 
& White 2001; Yoshida, Sheth, & Diaferio 2001), a semi- 
analytic model of galaxy formation (to track gas cooling, 
star formation, and feedback from supernovae; Kauffmann 
et al. 1999; Springel et al. 2001) and the Monte Carlo ra- 
diative transfer code CRASH (to follow the propagation of 
ionizing photons into the IGM; Ciardi et al. 2001; Maselli, 
Ferrara, & Ciardi 2003). The linear box size is 20h-^ Mpc 
comoving, which corresponds to 10.59 arcmin at 2: = 9.26. 

Simulations were run with two sets of parameters for 
the stellar emission properties. The run labeled S5 {late 
reionization model) in CFW adopts a population of PopIII 
stars with a Salpeter Initial Mass Function (IMF) and an es- 
cape fraction of ionizing photons fesc = 5%. The L20 {early 
reionization) run adopts a mildly top-heavy Larson IMF and 
fesc = 20%. Complete reionization (i.e. am ionized fraction 
per unit volume > 0.99) is reached by a w 8 and z « 13, 
respectively. The output of the simulation provides, among 
other quantities, a series of 128^ arrays with gas number 

* The simulations assume a ACDM "concordance" cosmology 
with a^=0.3, f2A=0.7, h =0.7, 05=0.04, n=l and <T8=0.9. 
These parameters are within the WMAP experimental error bars 
(Spergel et al. 2003). 
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density and ionization fraction. These were used by CM to 
derive (from the equations described in the previous section) 
the expected 21 cm emission from the IGM at different red- 
shifts. 



4 BASIC RADIOASTRONOMY EQUATIONS 

In this Section we will briefly introduce the concepts of radio 
astronomy relevant for this study (for more radio astronomy 
details see Taylor, Carilli, Perley 1998). The quantity mea- 
sured by a radio telescope is the sampled visibility function: 



K(w, v)S{u 



m) e:x.p[—2m(ul+i'm)]dl dm, (5) 



where {u, v) are the coordinates in the plane perpendicular 

to the viewing direction (Fourier or («, r) plane), v is the 
frequency of the observed monochromatic component, / (m) 
is the director cosine measured with respect to u (v) defining 
the object position in the sky, and S{u, v) is the sampling 
function (also called [u, v) coverage), defined by the projec- 
tion of the baselines vectors between the telescopes, which 
accounts for the fact that it is not possible to sample the 
visibility function on the whole {u, v) plane. The Fourier 
Transform of the sampled visibility function gives the dirty 
image: 



7^(Z, m) 



(w, v) S{u, v) exp[27ri(wZ -|- i'm)]dudv. (6) 



The dirty image is equivalent to the actual distribution of 
the radiation intensity, 1^, convolved with the dirty beam 
or point spread function (PSF). V,^ is measured in units of 
flux density (Jansky, Jy), while Ii, is a flux density per unit 
of solid angle (e.g. Jy/beam area) that is often described in 
terms of brightness temperature (K) through the conversion 
Tb = c^h/2k V. 

In a radio interferometer the monochromatic [u, v) cov- 
erage is given by the projections on the (w, v) plane of the 
baselines between the observing stations, such that 



sin HA^ 



cos 11 



- sin 5o cos HAq, sin (5o sin HAq 



(7) 



(8) 



where HAa and arc the hour angle and declination of the 
phase center, A is the mean band wavelength and L is the 
length of the baselines. Note that this approximation ignores 
any non-planar effects appropriate for small field imaging. 



5 LOFAR VIEW OF REIONIZATION 

LOFAR will explore the Universe at frequencies between 30 
and 80 MHz (Low Band) and between 110 and 240 MHz 
(High Band); it will be composed of thousands of anten- 
nas divided initially in 77 stations spread over a ^ 100 km 
radius, and eventually in up to 100 stations, with a max- 
imum baseline of ~ 400 km. This could be later extended 
to more than 1000 km if, in addition to the Netherlands 



where LOFAR is presently under construction, other Euro- 
pean countries join the project.^ The ~ 400 km baselines will 
guarantee a spatial resolution as high as « 1 arcsec at 150 
MHz. Several problems associated with low radio frequency 
observations need to be addressed, such as the ionospheric 
scintillation, interferences from human-generated radio sig- 
nals and extended/point source foregrounds. The sensitivity 
will be optimized for the reionization experiment by concen- 
trating about 40% of the collecting area in the inner (~ 2 km 
diameter) radius or Virtual Core of the instrument. The LO- 
FAR correlator is such that it may be configured in several 
different modes to achieve differing numbers of beams on the 
sky and spectral resolutions. More specifically, the product 
between the number of channels, A^'^h, and the number of 
polarizations is constant and equal to the total number of 
digital signal paths. In addition, the total available band- 
width of 32 MHz can be divided typically into up to 8 dif- 
ferent pieces to give different beams on the sky. Thus trade 
offs can be made between these parameters to optimize the 
array for the desired observations, i.e. bandwidth may be 
traded for an increased number of beams or polarization 
can be traded for an increased number of spectral channels. 

From the LOFAR provisory station configuration (100 
stations with bsiselines of up to ~ 360 km, J. Noordam, pri- 
vate communication) we have computed the (m, v) coverage, 
which is extremely complete as can be seen from Figs. 1-2. 
Fig. 1, in particular, shows the instantaneous (u, v) cover- 
age for a declination So ~ n/2 while the box in the lower 
right corner shows a zoom of the {u, v) coverage between -5 
and 5 Kilo-A; Fig. 2 reports the (w, v) coverage for a 8-hour 
continuous observing time at <5o = 7r/6. In this case the cov- 
erage is very dense thanks to the earth rotation synthesis 
and to the configuration of LOFAR's stations. 



5.1 Synthetic radio maps 

To produce synthetic radio maps of the reionization process 
as seen by LOFAR we proceed as follows. We denote by 
Ic — 20/128 Mpc the comoving thickness of a slice of 
the simulation box, and with: 



1.7 



Vo.iMHz/ V 10 y 



0.15 



-1/2 



Mpc, (9) 



the comoving length corresponding to an observational 
bandwidth Avobg. The number of slices that corresponds 
to a fixed Avobs is ric = Lc/lc- For a given Avobs, we extract 
from each simulation a 2D map of brightness temperature 
by averaging the values of 5Tb on a number of slices ric- 

The maps are first Fourier-transformed and then they 
are convolved with the (u, v) coverage, S{u^ v), given in 
Eqs. 7-8 and shown in Fig. 2. Such convolution is neces- 
sary to take into account the limitations due to LOFAR's 
sampling function. We assume a monochromatic (u, v) cov- 
erage, i.e. we neglect the width of the tracks in the (u, v) 
plane. This represents a simplification, however we find that 
for the present study this approximation has a small im- 
pact on the results in comparison to the sensitivity limits 



t See e.g. www.mpifr-bonn.mpg.de/staff/wsherwood/LOFAR/ 
white.paper.oct6.pdf 
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Figure 1. LOFAR instantaneous (u, v) coverage for a declination 
So = The units on the x and y axis are Kilo-A, where A is the 
mean band wavelength. The box on the lower right corner shows 
a zoom of the (u, v) coverage between -5 and 5 Kilo-A. 



200 L I I 



100 - 



I 



I I I — I I I I — I I I I — I I I I — I— 




01- 



100 - 



~ZQO r I I I I I'll ■ ■ I ■ 'ill I ■ I ■ ■ I ■ ■ I'll ■ ■ I ■ L 

-200 -100 100 

K-lambda 



200 



Figure 2. LOFAR {u,v) coverage for 5o = 7r/6 declination angle 
after 8 hours integration time, with the same units as in Fig.l. 



of the instrument, because of the exceptionally complete 
(w, v) coverage due to the log-spiral arms configuration of 
LOFAR. Next, an inverse Fourier transform has been per- 
formed to return to the image plane including the estimated 
instrumental sensitivity on different scales and frequencies 
as given at the official LOFAR website ^. The sensitivity 
of an interferometer is oc {Auobs^t)~^^'^', we assumed an 
observational bandwidth Avobs = 128 kHz and an observa- 



■t http://www.lofar.org/science 



tion time At = 1000 hours. Finally, we have convolved the 
resulting image with a gaussian beam of w 3 arcmin. This 
choice represents the best compromise between LOFAR sen- 
sitivity (which decreases for smaller beams) and the size of 
the simulation box, which is limited to ~ 11 arcmin. On this 
angular scale LOFAR will likely be able to detect the 21 cm 
signal thanks to its large collecting area in the Virtual Core. 
We note that this is mathematically analogous to applying 
a taper in the (u, v) plane to give greater weighting to data 
from the shorter baselines, which is entirely consistent with 
the most critical part of the data being collected by the Vir- 
tual Core. 

The results of this procedure are presented in Fig. 3 
for the S5 model, i.e. the late reionization scenario in which 
the EoR is located at a « 8, the most promising in terms 
of 21 cm line detection from the pre-reionization IGM. The 
maps represent the 21 cm emission brightness temperature 
in logarithmic scale, both directly from the numerical simu- 
lation (i.e. before the convolution with LOFAR instrumen- 
tal characteristics) and in the synthetic LOFAR images; the 
calculations were performed using the simulation outputs 
at redshifts z = 10.6, 9.89 and 9.26 (corresponding lo v = 
122, 130 and 138 MHz respectively). We find a region that 
should be detectable at a 3-ct level, showing that if reion- 
ization occurred relatively late LOFAR should detect the 
signal from the neutral IGM. The resolution and sensitivity 
of the instrument allow the identification not only of the 
largely neutral IGM, but also of the HII holes produced 
by the ionizing sources. The brightest emission peaks at 
z = 10.6 have STj, = 35 mK (the corresponding value in 
the original simulated map is ~ 0.2 K), whereas the holes 
are typically a factor of 5 below that - a dynamical range 
large enough to allow their robust identification as LOFAR's 
sensitivity will be ~ 10 mK at these scales, frequencies, ob- 
servational bandwidth and integration time. Even more im- 
portantly, the three synthetic maps demonstrate that the 
21 cm line tomography can provide a superb tool to track 
the redshift evolution of the reionization process (the accu- 
racy of this measurement could be noticeably enhanced by 
cross-correlating the 21 cm emission with CMB secondary 
anisotropies, see Salvaterra et al. 2005). 

We have repeated the same analysis halving (doubling) 
the observational bandwidth, i.e. Avoba = 64 (256) kHz; this 
corresponds to a sensitivity decrease (increase) by a fac- 
tor \/2 with respect the case analyzed above. The bright- 
est structures in the synthetic maps at z = 10.6 have 
5Tb = 34 mK for Avobs = 64 kHz with a dynamical range 
of 28 mK. The worse sensitivity (14 mK) obtained for this 
bandwith results only in a 2-a detection. For the largest 
Avob3 = 256 kHz, instead the peak 5Tb ~ 33 mK, and 
the dynamical range is 21 mK. Hence, given a sensitivity 
of w 7 mK, this would represent a 3-a detection in anal- 
ogy with the Auobs = 128 kHz case. Thus, we conclude that 
using the median observational bandwidth is preferable as 



3 We note that a total observing time of 1000 hours can only 
be obtained through a scries of observations taking place during 
a few months. There will be slight shifts in («, v) coverages be- 
tween the different observations, resulting in an ever better (u, v) 
coverage as sketched for a single 8 hour observations. 
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Figure 3. Logarithmic brightness temperature, logfiTj/A'], maps 
(Unear size fa 11 arcmin) of the 21 cm emission for the late reion- 
ization model S5 at redshifts z = 10.6, 9.89, 9.26 from top to bot- 
tom, respectively. Left panels: maps obtained directly from the 
simulation, i.e. before convolution with LOFAR characteristics; 
Right: LOFAR synthetic maps. 

it allows to have a fine-grained IGM tomography without 
missing too much signal. 

The situation is less promising in early reionization 
models as the L20 we have studied. In this case, the rapid 
disappearance of the neutral IGM component negatively af- 
fects the detection in two ways: (i) the instrumental sensitiv- 
ity decreases rapidly at low frequencies (i.e. towards higher 
redshift) ; (ii) the sky brightness temperature increases dras- 
tically. In this case, for a /S.Voba ~ 128 kHz (correspond- 
ing to a sensitivity of about 80 mK), we find a peak value 
STi, = 56 mK, a dynamic range of 26 mK. This corresponds 
to 0.32 a, implying an extremely difficult detection of the 
neutral IGM in an early reionization scenario. 

Thus we conclude that, although one cannot exclude 
that even early reionization scenarios may lead to a positive 
LOFAR detection, the most favorable case would be one in 
which the EoR took place relatively recently. If this is the 
case we have shown that in principle LOFAR will be able to 
see reionization in its full action. 



6 CONCLUSIONS 

By convolving simulations of cosmic reionization and 21 cm 
line emission with the provisional characteristics of the ra- 
dio telescope LOFAR, we have shown that LOFAR has in 



principle the sensitivity to detect the signal from neutral hy- 
drogen in the IGM, at least for a late reionization scenario 
(EoR at « 8) that appears to be more appropriate after 
the release of the third year WMAP data analysis results. If 
reionization occurred earlier, the experiment is much more 
challenging. Nevertheless, even an upper limit on the neu- 
tral fraction at different epochs would represent a corner- 
stone result which could also strongly constrain more com- 
plex/double reionization histories. Our study is an attempt 
at making specific predictions for LOFAR based on state-of- 
the-art simulations, and it is affected by a number of limi- 
tations and approximations. Although our simulation box is 
one of the largest used so far in numerical reionization stud- 
ies, the considered volume of (30 Mpc)^ is only marginally 
able to capture the global characteristics of cosmic reioniza- 
tion. Larger boxes could be used, but information on small 
scale structure (clumping) is easily smoothed out and lost, 
resulting in an extremely coarse description of the recombi- 
nation processes. It is only very recently (Kohler, Gnedin & 
Hamilton 2005, Iliev et al. 2006) that dedicated investiga- 
tions has started to tackle this very difficult problem. 

It is beyond the scope of this paper also to deal with the 
crucial problems of ionospheric scintillation and foreground 
contamination (i.e. Galactic free-free and synchrotron emis- 
sion, unresolved extra-galactic radio sources, free-free emis- 
sion from ionizing sources, synchrotron emission from cluster 
radio halos and relics); rather, we consider the idealized con- 
dition in which the reionization signal is the only one present 
in the sky. A number of studies have discussed these compli- 
cations in some detail (e.g. Shaver et al. 1999; Oh & Mack 
2003; Di Matteo, Ciardi & Miniati 2004) and we refer the 
reader to those papers for more information. We note how- 
ever, the success of the experiment will depend on the capa- 
bility of effectively removing the foregrounds and of correct- 
ing for ionospheric variations (e.g. Gnedin & Shaver 2004; 
Zaldarriaga, Furlanetto & Hernquist 2004; Santos, Cooray 
& Knox 2005; Morales, Bowman & Hewitt 2005). 
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